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Reconstruction of the primordial fluctuation spectrum from the flve-year WMAP data 
by the cosmic inversion method with band-power decorrelation analysis 
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The primordial curvature fluctuation spectrum is reconstructed by the cosmic inversion method 
using the five-year Wilkinson Microwave Anisotropy Probe data of the cosmic microwave background 
temperature anisotropy. We apply the covariance matrix analysis and decompose the reconstructed 
spectrum into statistically independent band-powers. The statistically significant deviation from a 
simple power-law spectrum suggested by the analysis of the first-year data is not found in the five- 
year data except possibly at one point near the border of the wavenumber domain where accurate 
reconstruction is possible. 



PACS numbers: 98.70.Vc, 95.30.-k, 98.80.Es 



I. INTRODUCTION 

Inflationary cosmology [l], [3, explains the origin of 
cosmic structure on various scales in an unified way. The 
expansion history during the inflationary stage of the 
early universe is recorded in the spectrum of seed struc- 
ture which can be revealed by modern cosmological ob- 
servations such as the Wilkinson Microwave Anisotropy 
Probe (WMAP) observation of the cosmic microwave 
background (CMB) I, i, i, 0, i i K [H [H, [11 • The 
WMAP mission evaluated every multipole moment of the 
CMB anisotropy spectrum in a wide range of scales which 
is potentially a record of the inflationary expansion his- 
tory with high time resolution. It is a feasible challenge 
to probe the highly time-resolved behavior of the inflaton 
field (s) which drives inflation. 

Since the initial release of the WMAP results [4, [sj, 
0) 0, Hi J it has been argued that the CMB tempera- 
ture anisotropy spectrum has nontrivial features such as 
running of the spectral index, oscillatory behaviors on 
intermediate scales, and lack of power on large scales 
[l3,[H,lilj[I3|i which cannot be explained by a power-law 
primordial spectrum that is a generic prediction of sim- 
plest inflation models. These features may provide clues 
to unnoticed physics of inflation. Some of these anoma- 
lous structures disappeared on the three-ye ar sp ectrum, 
however several anomalies are still existing jld, [lH. To 
explain these features, a number of inflation models have 
beenproposed H [H, [H, [H, H [H, [H, [13, [H, [H, 
[33,[3l|,|32,|3l[3j 

As an observational approach to these possible nontriv- 
ial features, which can be an alternative to model fitting, 
there have been several attempts to reconstruct the pri- 
mordial spectrum using CMB anisotropy data without 
any prior assumptions about the shape of the primordial 

spectrum [H, m [13, [H, [H, [10, [HI, [12, [H, O, [S [M [13 ■ 

One such attempt to reconstruct the primordial spectrum 
is a filtering method where the primordial spectrum is 
characterized by amplitudes on a few number of repre- 
sentative scales. While such methods have an advantage 



for reconstructing the global structure of the primordial 
spectrum, they may miss possible fine structures if their 
scale width is narrower than the filtering scale, which has 
been chosen rather arbitrarily so far. 

On the other hand, there exist other methods which 
can reconstruct the primordial spectrum as a continu- 
ous function without any ad hoc filtering scale to in- 
vestigate detailed features such as the cosmic inversion 
method [i^, [13, [U, [13l , the Richardson-Lucy method 
[53I . [5^ . [ssj , or a nonparametric method [s^ [53| ■ The 
cosmic inversion method has proved its ability of recon- 
structing the modulations off a power-law spectrum quite 
well by the analysis of mock anisotropy data. In the anal- 
ysis of the first-year WMAP data, we pointed out the 
possibility of nontrivial structures in the primordial spec- 
trum around the scales of 2.4 x lO^Mpc and 4.2 x lO^Mpc. 

Theoretically, according to the standard inflation 
paradigm, each wavenumber (fc-) mode of the power 
spectrum is mutually independent. On the other hand, 
each fc-mode of the power spectrum reconstructed from 
the observed CMB anisotropy has a strong correlation 
with neighboring modes because each multipole of CMB 
anisotropy, C^, depends on the fc-modes in the wide range 
around k = £/d where d is the distance to the last scatter- 
ing surface. In order to extract real features, therefore, 
it is important to decompose the reconstructed spectrum 
to mutually independent band-powers keeping resolution 
as high as possible. 

The purpose of this work is to apply the cosmic in- 
version method to the five-year WMAP temperature 
anisotropy spectrum [Tl| . update the reconstructed pri- 
mordial curvature spectrum, and perform the above- 
mentioned band-power analysis by diagonalizing the co- 
variance matrix. Then, we revisit the possibility of fine 
structures in the primordial spectrum. Because of the 
arbitrariness of the primordial spectrum, we inevitably 
incorporate infinite degree of freedom to our analysis, 
which results in degeneracy among spectral shape and 
cosmological parameters [ssl [59|. In this paper, we con- 
sider the concordance adiabatic ACDM model, where the 
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cosmological parameters (except for the ones character- 
izing the primordial spectrum) are those of the WMAP 
team's best-fit power-law model [13, , and instead fo- 
cus on the detailed sha pe of the primordial spectrum. 
Note that, as shown in [50], different choices of cosmo- 
logical parameters affect only the overall shape. The fine 
structures of the reconstructed power spectrum remain 
intact in the relatively small wavenumber range we probe. 

This paper is organized as follows: In Sec. II, the 
overview of our analysis is described. In Sec. Ill, we 
show the reconstructed primordial power spectrum from 
the five-year WMAP data and discuss its implication. 
Finally, Sec. IV is devoted to the conclusion. 



II. INVERSION METHOD 

A. Basic formulas 

Before presenting the inversion method we first list ba- 
sic formulas to be inverted. Although we only express for- 
mulas related to temperature anisotropy here, the same 
procedure can be repeated to polarization anisotropy as 
well [sTl , which will give us additional information in the 
future. 

The temperature anisotropy of photons coming from 
direction h observed at x is decomposed to Fourier modes 
and multipole moments as 
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where /i = k ■ n / \k\ = k ■ n and rj is the conformal 
time. In terms of multipole moment in the Fourier space, 
©^fcl??), which is defined by 
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aim is expressed as 



air,^iri,x)^ I ^e.(,7,M)^r;„(fc)e^^-. (3) 

Thanks to the assumption that there exist only adiabatic 
fluctuations, we can define the transfer function, Xi{k), 
from the Fourier mode of primordial comoving curvature 
perturbation, TZk{0), to Qekivo) by Oik = Xi{k)Tik{Q), 
where ?7o is the present conformal time. Then the an- 
gular power spectrum, Ci = (|a£mP), of temperature 
anisotropy and the power spectrum of curvature pertur- 
bation, P{k) = (|7^fc(0)|2}, are related by 
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This is the master equation we wish to invert. Note that 
Xi{k) depends on the cosmological parameters, too. 



B. Cosmic inversion method 

Let us introduce the cosmic inversion formula which 
relates the observational CMB anisotropy spectrum to 
the primordial curvature fluctuation spectrum by a first- 
order differential equation. Working in the longitudinal 
gauge, 

ds^ ^ a^irj) [-(1 + 2-i'{x))dr]^ + (1 + 2<P{x))dx^] , (5) 

the Boltzmann equation for 8fc(?7, fi) can be transformed 
into the following integral form [60j. 
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where the overdot denotes the derivative with respect to 
the conformal time. Here, and are the gauge- 
invariant quantities representing the Fourier transform 
of the Newton potential and the spatial curvature per- 
turbation in this gauge, respectively [6l|, [l^l, and 
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are the visibility function and the optical depth for 
Thomson scattering, respectively. In the limit that the 
thickness of the last scattering surface (LSS) is negligible, 
we find V{7'i) « <5(?7 — 77*) and e"'^^''-' « 9{ri — rj^,) where 
ry* is the recombination time when the visibility function 
is maximal [63|- Taking the thickness of the LSS into 
account, we have a better approximation for Eq. ^ as 

/•'(.end , 

where d = rjo ~ rj^ is the conformal distance from the 
present to the LSS and ?7*start and ?7*ond are the time when 
the recombination starts and ends, respectively. Here, we 
introduce the transfer functions, f{k) and g{k), defined 
by 
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We can calculate f{k) and g{k) numerically, which de- 
pend only on the cosmological parameters, for we are 
assuming that only adiabatic fluctuations are present. 
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Then, we find the approximated multipole moments as 

= UikMkd) + g{k)j',{kd)] n^io), (11) 

and the approximated angular correlation function as 



^crp,(i-^). (12) 



where C^^^ is obtained by putting Eq. into Eq. (g]), 
r is defined as r = 2dsin(6'/2) on the LSS, and inun and 
^max are lower and upper bounds on £ due to the limita- 
tion of the observation. In the small scale limit r <C rf, 
using the Fourier sine formula, we obtain a first-order 
differential equation for the primordial spectrum of the 
curvature perturbation, 

- k^f{k)P'{k) + [~2k^f{k)f{k) + kg^k)] P{k) 

= 47r / dr-^{r3C"^PP(r)}sinfcr = S(fc). 
Jo r or 

(13) 

Since f{k) and g{k) are oscillatory functions around zero, 
we can find values of P{k) at the zero-points of f{k) as 



P{ks) = 



ks g^iks) 



for f{ks) = 0, 



(14) 



assuming that P'{k) is finite at the singularities, k = kg. 
If the cosmological parameters and the angular power 
spectrum are given, we can solve Eq. ()13p as a boundary 
value problem between singularities. 

However, because Eq. ((T^ is derived by adopting the 
approximation ([8]), C'^'^^ is different from the exact an- 
gular spectrum for the same initial spectrum. The 
errors caused by the approximation, or the relative differ- 
ences between C^"^^ and Cf^ are as large as about 30%. 
Thus, we should not use the observed power spectrum 
Cf"" directly in Eq. Instead, we must find Cf^ 

that would be obtained for the real P{k). Although this 
is impossible in the rigorous sense, we have found an em- 
pirical remedy to find C'^^^ corresponding to Cf^^ in the 
following way. The crucial observation is that the ratio. 
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is found to be almost independent of P{k) |4§| . Using this 
fact, we first calculate the ratio, hf^ = C'^^'^^'^ /Cf^'^''^\ 
for a known fiducial initial spectrum P'"^ (fc) such as the 
WMAP team's best-fit power-law spectrum. Then, in- 
serting Cg^'^/bf^ , which is much closer to the actual Cj'^^, 
into the source term of Eq. ([13]), we may solve for P{k) 
with good accuracy. We may continue this procedure 
iteratively. 



C. Numerical calculation 

Given an initial condition and cosmological parame- 
ters, we can calculate the transfer functions /(fc) and g{k) 
numerically. Then, with the angular correlation function 
(or equivalently anisotropy spectrum), Ea. p3)) is solved 
as a boundary value problem between the neighboring 
singularities, and hence P{k) is reconstructed. Hereafter, 
we treat A{k) = k^P{k) instead of P{k) itself for the 
comprehensive display purpose and consistency with the 
common normalization of fluctuation amplitude. 

We adopt the adiabatic initial condition and fidu- 
cial cosmological parameter set which is the WMAP 
team's best-fit power-law model [H, E^l to calculate 
the transfer functions. For the reconstruction from the 
five-year WMAP data, the cosmological parameters are 
h = 0.724, VLb = 0.0432, f^A = 0.751, = 0.249, and 
r = 0.089. In this case, the positions of the singular- 
ities given by Eq. are kd ~ 70,430,690,..., where 
d ~ 1.42 X lO'^Mpc. Around the singularities, the recon- 
structed spectrum has large numerical errors that are am- 
plified by the observational errors. Using various model 
power spectra to investigate the accuracy of our inversion 
formula, we found that we can achieve the reconstruction 
with good accuracy in the limited range 120 < /cd < 380 
or 8.5 X lO-^Mpc"^ < A: < 2.7 x lO'^Mpc^^ In this re- 
gion, the errors due to our inversion method turn out to 
be much smaller than those due to observational errors 
including the cosmic variance. 

In practice, we cannot take the upper bound of the 
integration in the right-hand side of Eq. to be in- 
finite. The integrand in Eq. (fT3|) is oscillating with its 
amplitude increasing with r. To evaluate the right-hand 
side of Eg. p^ as finite, we convolve an exponentially 
decreasing function with a cutoff scale r^ut into the inte- 
gration of the Fourier sine transform. As the cutoff scale 
is made larger, the rapid oscillations of the integrand with 
increasing amplitude become numerically uncontrollable. 
On the other hand, if the cutoff scale is made smaller, the 
resolution in the fc-space becomes worse as Afc ~ 7r/rcut- 
For both numerical stability and resolution in fc-space, 
we adopt the optimized cutoff scale of rcut — 0.8d. 

For implementing the inversion scheme, we employ the 
routines of CMBFAST code 63] to calculate the transfer 
functions, but we have modified it to adopt much finer 
resolution than the original one in both k and £ so that we 
can compute the fine structure of angular power spectra 
accurately. 



D. Monte-Carlo simulation 

In order to incorporate observational errors and the 
cosmic variance and to obtain mutually independent 
band-powers, we employ Monte-Carlo method to calcu- 
late the covariance matrix of the reconstructed power 
spectrum. Producing 50000 realizations of a tempera- 
ture anisotropy spectrum based on the WMAP team's 
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best-fit power-law model whose statistics obey the likeli- 
hood function provided by WMAP team with good pre- 
cision, we obtain 50000 realizations of a reconstructed 
primordial spectrum. (The prescription of simulating 
anisotropy spectra is described in Appendix.) About 
1000 samples are sufficient for convergence of the covari- 
ance matrix introduced below, which assures the robust- 
ness of our conclusion. 

The covariance matrix of the reconstructed power spec- 
trum is defined by 



whose correlation matrix is diagonal and reads 



I ^ 1 ^ 1 ^ 

r? ^ Aa{ki)Aa{kj) — — -: ^ Aa{ki)— ^ ^^{kj) 



13=1 



^ {(A^ikMk,)))^ - ((A„(fc.)L {{A0{k,)})p , (16) 

where Aa{ki) represents the value of the reconstructed 
power spectrum at k = ki in the a-th realization, and 
Af = 50000 as mentioned above. The resultant Kij is 
a square matrix with its dimension equal to the number 
of sampling points ki, N. When we solve the cosmic 
inversion equation, (fT3|) , we discretize the relevant range 
of k to more than 2400 points. We calculate Ku at each 
point to estimate the error of A{k) there. In practice, 
however, the neighboring fc-modes are strongly correlated 
with each other, as mentioned above, and so the number 
of independent modes are much smaller. Hence we do 
not need to, and in fact, should not take so many points 
to calculate Kij. 

Since X is a real symmetric matrix, it can be diago- 
nalized by a real unitary matrix U, to yield 



C/X[/t = diag(Ai,A2,...,Ajv)=A, 



(17) 



where A^'s are the eigenvalues of K. We find they are pos- 
itive definite as they should be, provided that we take N 
small enough so that neighboring modes are not degen- 
erate with each other. In the present case, we find that 
if we take iV < 50, the covariance matrix is well behaved 
in the sense that the following procedure is possible with 
positive definite Xi and the well-behaved window matrix 
defined below. In terms of 



A- 



-V2^diag(A^^/^A-l/^...,A-^/^), (18) 



we define inverse square root of K as K = 
[/tA-i/2[/. 

We also define a window matrix W by 



(19) 



which satisfies the normalization condition J2j=i ^ij — 
1. Convolving Paiki) with this window function, we de- 
fine a new statistical variable Sa(ki) as 
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where *W denotes a transposed matrix. 

Note that in the previous band-power analysis of the 
power spectrum in the literature, decomposition into 
band-powers or wavelets is done by hand without cal- 
culating the covariance matrix; hence they result in ei- 
ther undersampling to lack traceability of fine structures 
or oversampling which generates unwanted correlations 
between different modes. 



III. RESULTS OF RECONSTRUCTION 
A. Reconstructed primordial spectrum 

Figure [1] shows the primordial spectrum A(k) recon- 
structed from the five-year WMAP data by the cosmic 
inversion method. In this figure, the solid wavy curve 
depicts the result of reconstruction and the solid straight 
line is the best-fit power-law ACDM model obtained by 
the five-year WMAP observations, namely, the power- 
law spectrum with A = 2.39 x 10~^ and Ug = 0.961, 
where A is the amplitude of curvature perturbation at 
ko — 0.002Mpc^^. The dotted curves around the power- 
law are associated la standard errors which correspond 
to the diagonal elements of the error covariance matrix 
(Fig. [3]) calculated by Monte Carlo simulation. 

The modulations of the reconstructed spectrum 
roughly fit inside the Icr borders; therefore this oscilla- 
tory deviation from the power-law spectrum is attributed 
to the fact that we are analyzing only one random real- 
ization of an ensemble of a simple power-law spectrum. 
It is not necessarily required that the inflation model re- 
sponsible for creation of our Universe should predict such 
a highly nontrivial spectrum, for our Universe is merely 
one of the realizations of quantum ensemble that is ac- 
companied by significant fluctuations. 

Once the possibility of prominent dips around kd ~ 
180, 220 and 350 was pointed out in the analysis of first- 
year WMAP data M- As seen in Fig. [H the modu- 
lations on such scales are fairly degraded in the recon- 
structed spectrum from the five-year data. Since many 
of the glitches and bites seen in the first-year WMAP 
anisotropy spectrum have disappeared in the five-year 
WMAP spectrum, the updated primordial spectrum is 
much more smoothed. 



B. Restored anisotropy spectrum 

Figure [2] illustrates the CMB temperature anisotropy 
spectrum which we restored from the updated primor- 
dial spectrum by adopting the cosmological parameters 
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FIG. 1: The reconstruction of the primordial spectrum by 
the cosmic inversion method from the five-year WMAP data. 
The sohd wavy curve is the result of reconstruction and the 
straight line is the best-fit power-law spectrum. Dotted lines 
are the standard deviation around the best-fit power-law. 
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FIG. 2: The comparison of the binned five-year WMAP data 
of = 4 with the anisotropy spectrum restored from the 
reconstructed primordial spectrum. 



of the best-fit power-law model and grafting the best- 
fit power-law spectrum outside of the investigated range. 
The effective value in the range of 120 < £ < 380 
for this restored anisotropy spectrum, which we calculate 
using the likelihood tool provided by the WMAP team 
[64], is 246, while that for the best-fit power-law model is 
273. Although the degree of fit is improved significantly, 
it does not have the original statistical meaning because 
we have incorporated a functional degree of freedom. 



C. Band-power analysis 

Despite the fact that our reconstructed power spec- 
trum restores the fine structures of the observed angular 



power spectrum of CMB anisotropy well, the oscillations 
observed in the raw reconstructed spectrum may not have 
true statistical meaning. As mentioned in the previous 
section and we can see in the covariance matrix Fig.[3]ex- 
plicitly, neighboring fc-modes of the reconstructed power 
spectrum are correlated with each other up to the range 
Afcd ^ 10. The origin of the correlation is the convolu- 
tion with the transfer function Xi(k) between A(k) and 
Cf. In particular, mapping the primordial spectrum to 
Ci erases the information that is responsible for the fine 
structure whose characteristic scale is much below the 
correlation width of Afcd ~ 10. 

It is important to decompose the power spectrum into 
mutually independent modes for evaluating the statisti- 
cal significance of the oscillations in the reconstructed 
spectrum. Here we construct band-powers using the 
window matrix W defined in Sec. IIIDl which diago- 
nalizes the covariance matrix and gives mutually inde- 
pendent errors. We take the dimension of the window 
matrix to = 40 and decompose A{k) in the region 
100 < kd < 400 to 40 band-powers. Figure S] shows the 
window functions for each mode. 

Figure [5] is the result of band-power analysis of the 
five-year WMAP data. In this graph, i-th data point 
indicates the value of 
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and the vertical error bar represents the variance 
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Here ki is the location of the peak of the i-th line of the 
window matrix Wij. The horizontal bar, on the other 
hand, indicates the width of the window matrix, where 
the dispersion of the fitted Gaussian is shown (see Fig. 
S]). We find their typical full width is Akd w 10 and 
the neighboring horizontal bars barely overlap with each 
other. As is seen clearly in Fig. [5l our band-power recon- 
struction of the five-year power spectrum basically agrees 
with the best-fit power-law as a whole. 

The i-th band-power depends on the multipole mo- 
ment on the relevant angular scale of £ ~ kid and also 
on the surrounding multipoles of £ ^ kid± 5. Indeed, we 
found that the statistics of each band-power subject to 
Gaussian distribution due to the superposition of several 
multipole moments even though the distribution of Ce 
is non-Gaussian (see Appendix). By virtue of our band- 
power analysis, we can also estimate the statistical signif- 
icance of the reconstructed spectrum itself by evaluating 
the deviation from the best-fit power-law spectrum at 
every band simultaneously. We estimated the whole sta- 
tistical significance of the central 34 bands which corre- 
spond to the scales of 120 ^ kd < 380 and found that, in 
terms of reduced value which is 0.64, the reconstructed 
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powers. The dotted line is a Gaussian fitting of tlie i = 10 
mode. 



FIG. 3: Tfie error covariance matrix of the reconstructed pri- 
mordial spectrum, 10^ A{k), calculated by Monte Carlo sim- 
ulation based on the cosmic inversion method. The central 
strip is indicating positive correlation and the surrounding 
strip is indicating negative correlation. Error enhancement 
on the largest scales of our analysis comes from the numeri- 
cal instability due to the singularity of transfer function while 
error regain on the small scale end is due to the domination 
of measurement error. 

spectrum from the five-year WMAP data is (— )1.5cr real- 
ization, ft is quite consistent with the best-fit power-law 
spectrum. 

Here the absence of the large modulations around kd ~ 
180, 220 and 350 is confirmed again. In the band-power 
analysis, we find dips around kd ~ 180, 220 and 350 
with the statistical significance about la. The nontrivial 
features reported previously have practically disappeared 
in the reconstructed spectrum from the five-year data. 

Note, however, that in Fig. [5l we find a prominent 
deviation around feci c^i 120 which would be a true signal 
of deviation from a simple power-law spectrum if our 
reconstruction method could be trusted there. While it 
is possible that this peak is indeed a real signal associated 
with the feature observed around £ ~ 120 in the angular 
power spectrum of the temperature anisotropy, it may 
simply be an artifact because fed ~ 120 is so close to the 
singularity fc^ ~ 70 that we can count on our method 
only marginally there. 



IV. CONCLUSION 

Using the cosmic inversion method, we have recon- 
structed the primordial power spectrum of curvature per- 
turbations assuming the absence of isocurvature modes 
and the best-fit values of cosmological parameters for the 
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FIG. 5: The band-power decomposition of the primordial 
spectrum obtained by the cosmic inversion method from the 
five-year WMAP data. The horizontal bar indicates the ef- 
fective width of each window function which depicts the dis- 
persion of the fitted Gaussian distribution. 

power-law ACDM model. While the range of accurate re- 
construction is rather narrow, 120 < fed < 380, we can 
reproduce the fine structures of modulations off a simple 
power-law with which we can recover the highly oscilla- 
tory features observed in C^. 

The statistical significance of the oscillatory structures 
in the reconstructed power spectrum is difficult to quan- 
tify due to the strong correlation among the neighbor- 
ing fc-modes. We have therefore performed the covari- 
ance matrix analysis to calculate the window matrix, 
W, which diagonalizes the covariance matrix into sta- 
tistically independent modes. We have chosen the large 
enough number of band-powers, N = 40, to probe the 
feature of the primordial spectrum as precisely as pos- 
sible while keeping the overlap of the window functions 
Wij small enough. 
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As a result we have found all the independent modes 
are consistent with the best-fit power-law ACDM model 
for the five-year WMAP data except possibly for a point 
around kd ~ 120. Whereas, the possible deviation from 
a simple power-law around kd ~ 180, 220 and 350 re- 
ported in the analysis of the first-year WMAP data [H^l 
has disappeared. This difference is due to the fact that 
the observed Cg in the five-year WMAP data has become 
much smoother than the first-year counterpart around 
the first peak. 

On the other hand, there still remain some nontriv- 
ial deviations from the expected Ci of a simple power- 
law spectrum around £ w fed w 40 even in the five-year 
data. Unfortunately, the cosmic inversion method can- 
not probe the primordial spectrum in this region. So we 
need to develop a different method which can probe the 
power spectrum for smaller wavenumber region kd < 120 
precisely. 

From the covariance matrix analysis, we have shown 
that statistically independent bands of the primordial 
spectrum have an effective width Akd « 10. This means 
that it is difficult to probe the possible fine structures on 
A{k) predicted by, say, trans-Planckian processes [65| . 
if the scale width of their characteristic modulation is 
narrower than Akd w 10. In other words, even if our re- 
sult is consistent with a simple power-law spectrum, this 
docs not rule out such possibility of narrow modulation 
of Akd < 10. 
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APPENDIX: PRESCRIPTION OF GENERATING 
MOCK ANISOTROPY SPECTRA 

In order to perform a proper error analysis of the re- 
constructed power spectrum, we must prepare a number 
of realizations of C^'s which obey the correct statistical 
distribution function. In the ideal situation with full-sky, 
uncontaminated observation, each multipole coefficient 
a£m is Gaussian distributed and mutually independent 
if primordial perturbation is random Gaussian. Conse- 
quently each angular multipole Ci is x^— distributed with 
2^-1-1 degrees of freedom, and different multipoles are 
uncorrelated. In practice, however, reliable observation 
can be made in only a finite fraction, fsky, of the sky 
and different £— modes are somewhat correlated. In this 
situation, the following form of likelihood function has 
been proposed [9|] and used in the statistical analysis of 



WMAP data 

-21n/:=i {Cf -Cf)Fu,{C? ^Cl'?) 

+ \ 5 - Zf^)F'u,{zF - Zl^l (A.l) 

where Cg = £{£ + \)CiI2t: and Zg = ln(C£ Mi) with 
Ml = i{£ + \)Ni/2iT. Here Fit is the Fisher matrix of Ci 
and Ff^i is related with Fu' as 

F'u, = (Cf + Mi)Fii, [Cf -f Ml,). (A.2) 

Superscripts D and th denote data and the theoretical 
model, respectively. 

The value of each component of Fisher matrix Fit 
is given in [(Mil where we find diagonal components are 
larger than neighboring off-diagonal elements typically by 
a factor of ^ 10^ — 10'^. When we depict C^'s, these diago- 
nal components of the Fisher matrix are used to indicate 
error bars associated with the respective multipole. How- 
ever, since different multipoles are correlated we would 
obtain an erroneous result if we created random samples 
based on these errors only. We should first find a basis 
consisting of linear combinations of C^'s which diagonal- 
ize the Fisher matrix to obtain statistically independent 
quantities. 

Here we describe the prescription to find an appropri- 
ate basis to diagonalize the Fisher matrix. First we note 
that, since the Fisher matrix is a real symmetric matrix, 
it can be diagonalized by a real unitary matrix which we 
denote by Q, namely, 

QtFQ = diag(r,„,,„,r,„,,„+i, • • • ,r,_J = i?, (A.3) 

where r^'s are eigenvalues of Fit and they are all positive. 
In terms of 

i?V2,diag(r;f^,r;/^^„...,r;fj, (A.4) 
we define matrices 

V EE QR^I^Q^ = V\ (A.5) 

and 

Vd = diag('i;£,„j„<>„^j„, • ■ • , (A. 6) 

that is, y is a "square root" of the Fisher matrix with 
= F, and Vd is a diagonal matrix with their diagonal 
components identical to those of V. We find 

{V-^Vd)^F{V-'Vd) = Vl 

that is, F can also be diagonalized if it is sandwiched by 
(y^iVo)t and V~^Vd- Note, however, that v^^'s are not 
eigenvalues of F because V~^Vd is not a unitary matrix. 
Nevertheless, defining 
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(A., 



and 



_ yth \ 



) / 



(A.9) 



we find the likelihood function is diagonalized as follows. 



) • (A. 10) 



Now Bi and are independent of each other. On 
the other hand, BiS are dependent on BiS. Reflect- 
ing the properties of the Fisher matrix and thanks to 
the normalization by Vj^^, we find that each diago- 
nal component of V^^V is equal to unity. Among the 
off-diagonal components in each line, {V^^V)u±2 have 
the largest magnitudes and their magnitude is no larger 
than 0.01. The other off-diagonal components are even 
smaller. Therefore, with a good approximation, we may 
put Be ^ Cf - Cl'' and Bi ^ {Cf' + J\fi)iZ^ - Zf"). In 
fact, we do not adopt these approximations but assume 
a relation 



B, 



C 



th 



(A.ll) 



which is inferred from the above approximate relations, 
so that each multipole is completely isolated in the like- 
lihood function (jA.lOp . With this approximation we can 
generate random numbers for each £ separately obeying 
a probability distribution function (PDF) 



'P[Bt\d\nBe oc exp 



(Cf +AA,)2ln2 1 + 



Bi 



dlnBt. (A. 12) 



Without the approximation (jA.lip . it would be compu- 
tationally formidable to generate many random samples 
of an angular power spectrum with the appropriate sta- 
tistical distribution. 

We generate random numbers for each B^ satisfying the 
above PDF (|A.12p from which we constitute realizations 
of Ci according to 





+1 




( 


+1 


( c\\ 

rth 

h 


n 




j 






„ ) 







(A.13) 



which have the desired correlation properties and PDFs. 

We make 50000 random realizations based on the 
WMAP team's best-fit power-law angular spectrum ac- 
cording to the above prescription and perform inversion 



using the cosmic inversion method to calculate the error 
covariancc matrix of the reconstructed primordial spec- 
trum. 
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